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Abstract. 

We study the scaling of the traces of the integer powers of the partially transposed reduced 
density matrix Tr(p^^ )" and of the entanglement negativity for two spin blocks as function of 
their length and separation in the critical Ising chain. For two adjacent blocks, we show that 
tensor network calculations agree with universal conformal field theory (CFT) predictions. In 
the case of two disjoint blocks the CFT predictions are recovered only after taking into account 
the finite size corrections induced by the finite length of the blocks. 
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1. Introduction 

During the last decade, it became clear that new concepts from quantum information theory 
can be extremely useful also to characterize many-body quantum systems both in and out of 
equilibrium. The most successful of these concepts is surely the entanglement entropy which 
turned out to be an optimal probe to distinguish the various states of matter, in particular 
for critical and topological phases (see e.g. [1] for reviews). Let p be the density matrix 
of an extended quantum system, which we take to be in a pure quantum state |^), so that 
p = |\E')(\1'|. Let the Hilbert space be written as a direct product T-L = Ha '^'Hb- ^'s reduced 
density matrix is = Tr^p and from this the Renyi entanglement entropy is defined as 

4") = ^lnTrp:^, (1) 

that for n = 1 reduces to the more studied von Neumann entropy. However, the knowledge 
of the Renyi entropies for any n provides more information than the n = 1 case and, in fact, 
it gives the full spectrum of the reduced density matrix [2] . 

For a one-dimensional critical system whose scaling limit is described by a conformal 
field theory (CFT), in the case when A is an interval of length i embedded in an infinite 
system, the asymptotic large i behavior of the Renyi entropies is [3, 4, 5, 6] 

S^:^ = Ul + -)log^- + c'^, (2) 
\ n / a 

where c is the central charge of the underlying CFT [7] and a the inverse of an ultraviolet 
cutoff (e.g. the lattice spacing). The additive constants are non universal but satisfy some 
universal relations [8]. 

Unfortunately, the entanglement entropies do not provide any information about the 
multipartite entanglement in a many-body system. In order to give a practical example, let us 
imagine to divide an extended quantum system in three spatial parts which we call Ai, A2 and 
B. We can define the reduced density matrix PaiuA2 = T^^bP, but the corresponding entropies 
^AiuA2 would measure only the entanglement between A = Ai U A2 and B, but not the 
entanglement between Ai and A2 which is the tripartite entanglement in the system. Even the 
Renyi mutual information Iai-A2 ~ + ~ '^a^iuA2 ^'^^^ ^'^^ provide the entanglement 
between Ai and A2, but gives only a measure of their correlations (see e.g. Ref. [9]). The 
quantification of the tripartite entanglement in a pure state has been longly a problem (see e.g. 
Refs. [1,9, 10, 11, 12]) until a computable measurement of entanglement has been introduced 
by Vidal and Werner [13]. This measure has been called negativity and is defined as follows. 
Let us denote by |e-^'') and \ef^) two bases in the Hilbert spaces corresponding to Ai and A2 
respectively. Let us define the partial transpose with respect to A2 degrees of freedom as 

and from this the logarithmic negativity as 

^ = ln||PAWJI=lnTr|pJu^J, (4) 

where the trace norm Hp^^^u^^H '•^^ absolute values of the eigenvalues Aj of 

PaiuA2- "^^^ negativity has been used to investigate the tripartite entanglement content of 
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many body quantum systems both in their ground-state [14, 15, 16, 17, 18, 19, 20, 21] or out 
of equilibrium [22, 23]. Furthermore, the definition of the negativity is very intriguing because 
it is basis independent and so calculable by quantum field theory (QFT). Only recently [24, 25] 
a systematic and generic method to calculate the negativity in QFT (and in particular CFT) 
has been developed. This method is based on the calculation of the traces Tr(p^^)" and the 
negativity is recovered in a replica limit. The negativity between two adjacent intervals has 
been calculated for a general CFT and turned out to be a universal quantity depending only on 
the central charge. The negativity between two disjoint intervals is still universal but depends 
on the full operator content of the theory and the corresponding moments Tr(p^^^^^^)" have 
been calculated explicitly for the free compactified boson [25]. In this manuscript we extend 
the previous results to the critical Ising CFT and confirm our predictions by explicit tensor 
network calculations in the transverse field Ising chain described by the Hamiltonian 



where crj'^ are Pauli matrices acting on the spin at site j and we use periodic boundary 
conditions. The model has a quantum critical point at /i = 1 and in the continuum limit 
corresponds to a free massless Majorana fermion which is a CFT with central charge c = 1/2. 

The manuscript is organized as follows. In Sec. 2 we discuss the general QFT approach 
to negativity [24, 25] and report those results which are valid for an arbitrary CFT. In Sec. 
3 we report the CFT calculation of the integer powers of the partial transpose of the reduced 
density matrix for two disjoint intervals in the critical Ising theory. Although the Ising chain is 
mapped to free fermions, it is still not known how to calculate effectively the partial transpose 
and for this reason we resort purely numerical methods. In Sec. 4 we introduce the tree tensor 
network approach and explain how to adapt it to the calculation of the eigenvalues of the 
partially transposed reduced density matrix. Finally we analyze the numerical data in Sec. 5 
taking properly into account the corrections to the scaling. In Sec. 6 we draw our conclusions. 

2. General results for the negativity in CFT: one interval and two adjacent intervals 

In the following we will be interested in tripartion of a one-dimensional systems such as those 
depicted in Fig. 1. We will denote with the reduced density matrix of A = AiU A2, i.e. 
Pa = which is obtained by tracing out the part B of the system, i.e. pa = Tr^p. 

The quantum field theory approach of negativity is based on a replica trick [24, 25], i.e. on 
the calculation of the traces Tr(p^^)" of integer powers of pj. According to Eq. (4) we are 
interested in the sum of the absolute values of the eigenvalues Aj of p^. It turns out that 
Tr(p^^)" have a different functional dependence on |Ai| according to the parity of n. Indeed, 
for n even and odd (that we denote as Ue and Uo respectively), the traces of integer powers of 

To 

are 



L 




(5) 




(6) 



Ai>0 



Ai<0 
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If we set Tie = 1 in Eq. (6) we formally obtain Tr|p^^| in which we are interested. Instead, 
setting no = 1 in Eq. (7) gives the normalization TrpJ = 1. This means that the analytic 
continuations from even and odd n are different and the trace norm in which we are interested 
is obtained by considering the analytic continuation of the even sequence at Ue — )■ 1, i.e. 

S = lim lnTr(p^2)"'= . (8) 

rie— S>1 

In a QFT, the traces of integer powers of the partial transpose for two disjoint intervals 
(as in Fig. 1 (top)) are partition functions on n-sheeted Riemann surfaces or equivalently the 
correlation functions of four twist-fields [25] 

Tr(pJ)" = {rn{u,)rn{Vl)Tn{u2)rn{v2)) , (9) 

i.e. the partial transposition has the net effect to exchange two twist operators compared to 

Tipl = {Tn{Ui)Tn{Vl)Tn{u2)Tn{v2)) • (10) 

(See Refs. [26, 6] for an introduction to the concept of twist fields.) 

Eq. (9) is of general validity, but it simplifies when specialized to the case of two adjacent 
intervals, obtained by letting vi — )■ M2, giving the three-point function 

Tr(pJ)" = {Tn{u,)Tl{u2mv2)) . (11) 

A further simplification occurs when specializing to a pure state by letting S — )• (i.e. 
U2 — )■ vi and V2 — i- Ml) for which Tr(p^^)" becomes a two-point function 

Tr(pJ)" = (r„^(«2)r'M). (12) 

As explained in more details in Ref. [25], as a partition function on an n-sheeted Riemann 
surface, this expression depends on the parity of n because connects the j-th sheet with 
the (j + 2)-th one. For n = rie even, the ng-sheeted Riemann surface decouples in two 
independent (ne/2)-sheeted surfaces. Conversely for n = Uo odd, the surface remains a no- 
sheeted Riemann surface. In formulas, these observations are 

Triplr- = {{%^/2{u2)Tn./2{v2))f = (Trp-/^)^ , (13) 
Tr(pJ)- = (J-nMr nA^2)) = TtpX ■ (14) 

Hence, for a bipartite system, Tr(p^^)" can be generically written as a function of Trp^^, as it 
should. In particular, taking the limit He — > 1, we obtain the logarithmic negativity £ = S^^^'^^ 
which is the well known result that for bipartite states the logarithmic negativity equals the 
Renyi entropy of order 1/2 [13]. 

2.1. A single interval in a CFT. 



For conformal invariant theories it is useful to consider first the case of a reduced density 
matrix corresponding to a bipartite system obtained by letting 5—7-0. We recall first the 
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standard result for Trp^^ in a bipartite pure state in the case when A2 is an interval of length 
i = |m2 — ^2! in the infinite line [4] 

/£\ -c/6(n-l/n) 

T^^PM = {rn{u2)Tn{v2)) = C„ f " j , (15) 

i.e. the twist fields behave like primary operators with scaling dimension A-j-^ = = 
c/12(n - 1/n). 

Following Refs. [24, 25], when the interval A2 is embedded in an infinite system, we can 
derive the powers of using Eqs. (13) and (14) and from the general formula for Trp^^ in 
Eq. (15), obtaining 

^ ,„ //\ -c/3(ne/2-2/ne) 

Tr(pJ)- = {{Tnj2{u2)rnj2{v2W = (T^pX^y = clJ-) , (16) 

\(J/ / 

and 

Tr(pJ)- = (r„„(«2)r„„(t;2)) = Tip- = c„„ (-) . (17) 

This simple result shows an important feature of the negativity in CFT, i.e. for n = Ue even, 
7? and T„ have dimensions 

Ar£=Ar^ =^f^--V (18) 



6 V 2 n. 



while for n = rio odd, 7? and T„ have dimensions 

^ 't-O 'I'D 

Ar2 =A^2 =-^fn„--y (19) 



■^no 12 V 

the same as 7^^. Thus, performing the analytic continuation from the even sequence, we 
finally have 

||pj||= limTr(p^^)"^ = c2/2(-)'^' £ = ^\n- + 2\ncy2, (20) 

which again is the result that for pure bipartite states the logarithmic negativity equals the 
Renyi entropy of order 1/2. Continuing instead the odd sequence to tt-o — >^ 1, we recover the 
normalization Trp^^ = 1. 

Notice that the constants Cn^_„^ are non-universal, but they are the same as in the 
entanglement entropies and so new non-universal constants do not appear in the partial 
transposition, but are all already encoded in the reduced density matrix. 



2.2. Two adjacent intervals in a CFT. 

Let us now consider the configuration of two intervals Ai and A2 of length £1 and £2 sharing 
a common boundary at the origin as graphically depicted in Fig. 1 (bottom). This can be 
obtained by letting t>i — )■ M2 = in Eq. (9) and it is then described by the 3-point function 

(we set ui = —ii and f 2 = £2) 

T^{PaT = {Tni-ilKmnm , (21) 
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Figure 1. The two configurations of the Ising chain that we consider. Top: the entanglement 
between two disjoint intervals Ai and A2 embedded in the ground-state of a larger chain. 



Bottom: The entanglement between two adjacent intervals in a larger system. 



whose form is determined by global conformal symmetry [27] 



where, for simplicity, we have set the UV cutoff a = 1. Notice that dn is not universal, but 
can be written as dn = clC^ ^2 where the structure constant ^2 is universal and can 
be determined by considering the proper limit of the four-point function as explicitly done in 
Ref. [25]. 

For n = He even, using the dimensions of the twist operators calculated above, we find 
the universal scahng relation 

that in the limit rig -> 1 gives 
For n = Uo odd we have 

T2 



that for no — )■ 1 gives Tip^ = 1 as it should 



2.3. Finite systems. 



All previous results are generalized to the case of a finite system of length L with periodic 
boundary conditions by using a conformal mapping from the cylinder (with axis perpendicular 
to the spatial coordinate) to the plane. The net effect of the mapping is to replace each length 
(.i with the chord length (L/n) sin^irii/L) in all above formulas. 
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Thus for the case of a pure state in a finite system the generalization of Eqs. (16) and 
(17) are 

2 (L 7^£^^-c/3K/2-2/ne) 



Tr(p^^)- 



Cne/2 I — Sm 



vvra L 

L 7r£\ -c/6(no-l/no) 

c„„ I — sin — 
^vra L / 

o = -ml — sm — 1 + 2 m ci/2 ■ 



Tr(pJ)"° 



(26) 
(27) 

_ , . (28) 

For two adjacent intervals embedded in a finite system, the generalizations of Eqs. (23) and 
(25) lead to the universal scaling relations 

Tiipl^r — 

Tr(pJ)- 



[g sin (^) sin, ^ 



^^-|c/6(ne/2-2/ne) 



TT L 



1 c/6(ne/2+l/ne) 



^4sm(^) sm(^) sin^(^ 



\^ L y "^^^^ V ^ / 
and for the logarithmic negativity 



c/12{no-l/no) ' 



(29) 



S = -ln 
4 



L sin (^) sin , ^ 



^ sin 



7r(^l+^2) 



+ cnst. 



(30) 



3. Negativity for two disjoint intervals in the Ising universality class 

For a generic CFT with central charge c, by using global conformal invariance, the four point 
function of twist fields can be written as [27, 28] 

c/6{ri— l/n) 



^31^42 



^21^43^41^32 

Zj and X is the four-point ratio 



(31) 



(32) 



{Tn{zi)Tn{z2)Tn{z^)Tn{Zi)) = C„ 

where the positions Zi G C, % = z, 

Z21Z43 
X = , 

2:31 2^42 

and X its complex conjugate. (All distances are measured in units of the UV cutoff a.) The 
function J-'nix, x) is real and must be computed case by case because it depends on the full 
operator content of the theory as explicitly done in some CFTs [28, 29, 30, 31, 32, 33, 34, 35] 
and also numerically in some lattice models [36, 37, 38, 39, 40, 41, 42], but mainly with x 
real and < x < 1. Since 2:3 1-242/ (-22 1-^43 -241^32) = l/(-22i-243(l — x)), the normalization of 
J-'nix, x) is given by J-'n(0,0) = 1. Thus, in the limit a; — )■ (i.e. large distance between 
the two intervals), the four-point function becomes the product of two two-point functions 
normalized through the constant c^. 

As explained in Refs. [24, 25], in order to give a physical interpretation to Eq. (31) as 
Trp^ or Tr(p^^)", the points Zi must be on the real axis and their order allows to distinguish 
between Trp^ and Tr(p^^)". We denote the boundaries of the two intervals with four real 
variables uj and vj (j = 1,2) such that ui < Vi < U2 < V2 in both cases. With these 
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definitions, Eq. (31) for Trp^ corresponds to zi = ui, Z2 = fi, = U2 and z^ = V2, leading 
to 

(/ \/ \ \ c/6(n— l/n) 

1 , ' ' ' 'I ,) -^«(-). P3) 

{Vi - Ui){V2 - U2){V2 - Ui){U2 - I'l) J 

where 

_ ^21^43 _ (^1 -'»l)(^2 -'»2) 
^31^42 [U2 - Ui)[V2 - Vi) 

and < s < 1. Being x real, the dependence on x = a; has been dropped. 

Instead, Tr(pJ)" is obtained from the four-point function (31) with the choice zi = ui, 
Z2 = Vi, Zs = V2 and Z4 = U2, namely by exchanging 2:3 ^ Z4 with respect to the previous 
case. The result can be then written as 

(/ N / \ \ c/6(n-l/n) 

? , ' 'I ;) SM, (35) 

with 

^21^34 (^^1 - Ui){v2 - U2) 

y = = 7 ^7 V , (36) 

^412^32 {U2 - Ui){V2 - Vi) 

and < y < 1 and again we dropped the dependence on y because y is real. Notice that we 
have chosen the definitions of x and y in such a way that the functional dependence on Uj and 
Vj is the same, but they are two different four-point ratios, as their dependence on Zi explicitly 
shows. The relation between the two is y = x/{x — 1). The expressions (33) and (35) are 
related as 

= (1 - y)^/3(- i/n) jr^ f_!L^\ , (37) 

\y - IJ 

as easily derived from Eq. (31). 

The logarithmic negativity is obtained by considering the even sequence n = and 
taking its analytic continuation Ue ^ I. Thus, the functional dependence on n of J-'nix) for 
X = y/iy — I) < must be different for even and odd n, in such a way that the limit 



£{y) = lim In Qn^{y) = lim In 



y 



(38) 



^y-i 

is not trivially equal to lim In Qudy) = 0. It is natural and convenient to introduce the ratio 



J.1 



R (y) ^ ^<PaT ^ {Tn{u,)Tn{Vl)Tn{u2)Tn{v2)) 



Gn{y) 



yy/3{n~l/n) ^nijj / {jj - I)) 



^n{y) J^niy) 

Indeed, in such ratio the prefactors cancel and a function only of y G (0, 1) is left. Moreover, 
the ratio in Eq. (39) is an easy quantity to consider in the numerical computations. Notice 
that, being lim Trp^ = 1, the logarithmic negativity is obtained equivalently from the replica 

n— >-l 

limit of this ratio 

S{y) = In lim R^M- (40) 

n.e->l 
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Thus, while Trp^ requires J-'n{x) only for x E (0,1), in order to extract Tr(p/)", and 
consequently the negativity £{y), one has to compute the function J-'nix) also for x < 0. 

We remark that for the special case n = 2, we have T2 = T2 because 7^ is the identity 
operator. This implies that i?2 (y) = 1 for < y < 1 identically. 

In order to find x) for the Ising model, it is useful to review the corresponding 

result for the free compactified boson (i.e. a Luttinger liquid field theory), derived in Ref. 
[25]. 



^) = , (41) 



3.1. Free compactified boson 

For the free real boson compactified on a circle of radius R, the function J-'nix, x) in Eq. (31) 
for any complex x has been computed in [25]. This has been possible starting from some 
results in Refs. [43, 44] and generalizing to x G C the procedure employed in Ref. [29] for 
< X < 1. The final resuk is [25] 

e(o|r,(x,x) ) 

where we introduced Fk/nix) = 2Fi{k/n, 1 — k/n; 1; x) to denote the special case of the 
hypergeometric function 2-^1 (a, b; c; z). Here 9(01-?^) is of the Riemann-Siegel theta function 
which is generically defined as 

e{0\K) = ^ exp[i7rm-ir-m], (42) 

me ZP 

where K is ap x p symmetric complex matrix with positive definite imaginary part and is 
the p dimensional vector made by zeros. The 2(n — 1) x 2{n — 1) symmetric complex matrix 
T^(x, x) entering in Eq. (41) is 

^"(^•^' - ( 'n rllr, ) • 

where the parameter i] is proportional to R"^ and the (n — l) x (n — l) symmetric real matrices 
TZ and I are respectively the real and the imaginary part of the (n — 1) x (n — 1) period matrix 

2 "^^ 

r(x) =TZ + iI= - ^ Tk/nix) sm{nk/n) Ck/n , (44) 



n 

k=l 



with 



_ . Fk/n{l - X) 

= ' F ' (45) 

^k/n[X) 

and Ck/n is the (n — 1) x (n — 1) symmetric matrix whose elements read 

(Cfc/n)„ = cos[27rA;/n(r - s)] , r, s = 1, . . . ,n - I . (46) 

Notice that, because of the k n — k invariance of the functions for n odd the 

denominator in (41) becomes the square of a product over k going from 1 to {uo — l)/2. 
Similarly, (44) can be written as a sum over k going from 1 to {uo — l)/2. This suggests that 
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the term corresponding to k/n = 1/2, occurring only for even n, plays a key role in the limit 
(38), as shown explicitly in [25] for the non-compactified boson in the hmit of close intervals. 

Few comments are now in order regarding Eq. (41) before proceeding to the partial 
transpose. Comparing Eq. (41) with old results about CFT on higher genus Riemann 
surfaces [44, 45, 46, 47, 48], one observes that J>j(x,x) is related to the partition function 
Z = Z'^"Z'^^{R) on a genus g = n — 1 Riemann surface, which can be factorized into a 
quantum part and a classical part, and all the dependence on i? of Z is contained in the latter. 
In particular, one finds that Z'^\R) = 0(O|T^(x,x)) and the period matrix of the Riemann 
surface is r(x) given in (44). As already remarked [43, 29], here we are not dealing with 
generic genus g Riemann surfaces, but with a subclass of them obtained through the replica 
method (see the appendix of [30] for a pictorial representation). Thus, while for Trp^ one 
needs the period matrix t{x) only for < a; < 1, where its real part TZ vanishes identically, 
for Tr(pJ)" the period matrix r(x) for x < is required and in this regime TZ enters in a 
crucial way. The expression (41) is explicitly invariant under t] ^ l/r]. We recall that one 
can also rewrite Eq. (41) in a form which makes manifest its invariance under x ^ 1 — x. 
This form is particularly useful also to study the regime — )■ oo of non-compactified boson 
because some simplifications occur and the Riemann- Siegel theta function does not contribute 
in this limit [25]. 

At this point, it is easy to write Tr (p^^ )" for the compactified boson from Eq. (35), where 
Qn{y) is given by (37) with c = 1 and J^n{x) by (41), i.e. 



where < y < 1. From Eqs. (47) and (41), we can write the ratio (39) for the free 
compactified boson 



In the decompactification regime r/ — )■ oo, this formula has been checked numerically in 
Ref. [25] for an harmonic chain with periodic boundary conditions, whose Hamiltonian is 
the lattice discretization of a free boson (Klein-Gordon action). Always in this regime and in 
the limit ?/ — )■ 1^ of close intervals, the analytic continuation rig — )• 1 of the even sequence 
Rneiy) has been explicitly carried out [25]. The analytic continuation rig — )■ 1 for finite and 
generic r] is still unknown. A similar problem occurs for the von Neumann entanglement 
entropy through the analytic continuation of the Renyi entropies of two disjoint intervals. We 
finally mention that the results for some finite i] have been recently checked in Monte Carlo 
simulations [50]. 

3.2. Ising conformal field theory 

In the previous section we have shown that J-'nix, x) can be found from the partition function 
of the corresponding model on the genus g = n — 1 Riemann surface obtained from the 
replica method. The partition function of the Ising model on a genus g Riemann surface has 



g^y) = {I - yY^-y-y^ 



(47) 




(48) 
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Figure 2. The function i?„ {y) for the Ising conformal field theory as a function of y for several 
values of n (n = 2, 3, 4, 5 curves from top to bottom). 



been derived in the eighties [45, 46, 48, 49] by considering Z2 orbifolds models on generic 
Riemann surfaces. From the knowledge of the partition function Z = Z'^^Z'^^(R) of the 
compactified boson on the same Riemann surface, we can explicitly access the period matrix 
on the Riemann surface. Then, using the mentioned old results for the Ising model on a genus 
g Riemann surface, we can finally write [30] 

Rising = iZ^n'^^2-aJ2memr)\ , (49) 



where r is the g x g period matrix. Here we introduced the Riemann- Siegel theta function 
with characteristic e, which is defined as 

gi7r[{m+e)-r-(m+£)+2(m+£)-(z+(5)] (50) 

where z G C^/ (Z^ + r Z^) is a g dimensional complex vector and the characteristic e is given 
by a pair of g dimensional vectors e and 6 made by O's and 1/2's 



e|el(z| 



E 



■ S ' 




- El,. 


■ ■ ,Bg 


. 6 . 






..,6g . 



e,,5,e {0,1/2}. 



(51) 



In the case at hand, the period matrix r(x) is given by Eq. (44) and the partition function for 
the compactified boson Z = Z'^^Z^\R) is Eq. (41). Thus, for the Ising model we have 

Ee|0[e](O|r(x))| 



rn-l 



lk=l \Fk/n{x)\^^^ 



(52) 



For < X < 1, this function reduces to J-'n{x) found in [30] for the Renyi entropies of the 
Ising model and it has been tested numerically through various methods [38, 39, 40], finding 
agreement once the finite size corrections are properly taken into account. 
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Writing Tr(p^^)" for the Ising model is now easy by using Eq. (35) with c = 1/2, and 
Qn{y) from Eqs. (37) and (52), namely 

where < y < 1 and r is given by Eq. (44). From this expression and Eq. (52), we can write 
the ratio (39) for the Ising model| 

R M - ,1 _ ,)(..-./«)/» fjlJkM. Ee|e|e|(0|r(^))| 

«„(!/) -(1 y) 1\\f„4^^) EeieN(o|rfe))| ' 

The curves Rn{y) for < y < 1 are shown in Fig. 2 for n = 2, 3, 4, 5. They are concave 
curves with -R„(0) = 1, which stay below 1 for n > 2. When n = 2, we have R2{y) = 1 
identically, as we will explicitly prove below. Moreover, Ruiiv) < RnAv) when n2 > rii, 
for all < y < 1. Approaching y = I, Rn{y) approaches zero very steeply in a calculable 
power-law way that, for a general model, depends only on the central charge [25]. 

It is worth mentioning that after a superficial analysis, one could erroneously conclude 
that the sum over the characteristics occurring in Eq. (52) involves 2^^ terms. This is not true 
because many of these terms are identically zero. In order to see this, let us introduce the 
parity of the characteristic e which is the parity of the integer number Ae ■ 6 which is related 
to the transformation property of the Riemann-Siegel theta function [46] 

e[e](-z|r) = (-l)^-^e[e](z|r). (55) 

There are 2^~^{2^ + 1) even characteristics and 2^~^{2^ — 1) odd ones. Eq. (55) implies 
that 0[e](O|r) = identically for odd characteristics. This means that the sum over the 
characteristics which defines J-'n(x, x) in (52) is made only by 2""^(2"^^ + 1) terms, namely 
the ones with even e. 



3.3. The special case n = 2 for Ising model. 

For the special case of n = 2, since T2 = T2, we have R2{y) = 1 identically in < y < 1. 
This is true for any CFT, as already noted above, and therefore also for the Ising model. 
Nevertheless, we find it instructive to see explicitly how this simplification occurs in Eq. (54), 
since this is not apparent from the formula and the explicit proof is related to the modular 
invariance of the torus partition function [51]. 

When n = 2 the Riemann surface is a torus and the period matrix reduces to its modulus 
r(x) = Tii2{x), which is a complex number. In this case the four-point ratio x G C and its 
modulus r(x) are related as [43, 28] 

F^/2{x) = -K{x) = e,{T{x))\ r(x) = ^^^^~^\ ^(^)=f^V' (56) 
TT K{x) V^3(r)/ 

X We have been informed that this resuh has been derived independently also by V. Alba [50]. 
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where K{x) is the complete elliptic integral of the first kind. From the properties of K, for 
< y < 1 one finds [52] 

1 \ 2 

TT 

2 

TT 



Fl/2 
Fl/2 



y 



y-i 



■^l^{K{l-y)-iK{y)), 



(57) 
(58) 



Given r = r(y) in the second formula of (56), (57) and (58), it is straightforward to observe 
that for < y < 1 we have 

y 



T = T 



y-i 



r + l. 



(59) 



We also need to recall the following transformation properties of the Jacobi theta functions 
under the modular transformation r — )• r + 1 [27, 51] 

92{t + 1) = e'^/%{T) , 03(r + 1) = 9,{t) , 9,{t + 1) = e^ir) . (60) 

From the first formula in (56), (59) and the second formula in (60), one finds that 



Fl/2 



y 



y-i 



9Jt) 



(61) 



Given the expressions above, we are ready to prove that R2{y) = 1 for < y < 1. Indeed, 
from Eqs. (56) and (61) we get 



y) 



1/4 







Fi/2{y) 


1/2 















(62) 



where in the first equation the identity 9f = 9^ — 02 has been employed. As for the ratio 
containing the Riemann-Siegel theta functions in (54), when n = 2 and < y < 1, it 
simplifies to 



\d2ir')\ + 




^4(r')| 


\02{r)\ + 


^3(r)| + 





1 



where Eq. (59) and the modular transformations (60) have been employed. 



(63) 



4. The partial transposition with Tree Tensor Networks 

The Ising chain in a transverse magnetic field is the most studied one dimensional model and 
this is mainly due to fact that it is a non-trivial model that can be mapped to a system of 
free fermions by means of a Jordan- Wigner transformation [53]. However, it is still a very 
hard problem to find an effective way to calculate the partial transpose of the reduced density 
matrix (or at least its eigenvalues) in the free fermions formulation as opposite to models of 
free bosons for which many results are already available (see e.g. [14, 16, 25]). For this 
reason, in order to check the analytical predictions from conformal field theory, we resort 
to purely numerical methods based on tree tensor techniques which already have been very 
effective in calculating of the entanglement entropies of two disjoint intervals [38, 40]. An 
alternative approach based on classical Monte Carlo simulations (generalizing the method in 
Refs. [54, 38, 40, 55]) has been also independently developed at the same time of this work 
by V. Alba [50]. 
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Figure 3. Graphical representations of some quantities for a pure state and of some useful 
properties in the language of the tensor networks. Top: We consider the bipartition of the 
pure state j^*) in Ai and A2 and its Schmidt decomposition through the isometrics Ui and U2- 
Bottom left: the density matrix p^. Bottom right: the partial transpose p^. 

4.1. Tensor networks: notation and examples 

We consider a ID lattice C made of L sites and on each site s G £ we set a local Hilbert 
space Ns of finite dimension d (e.g. for a spin-1/2 for which d = 2). An arbitrary pure state 
IM/) G V®^ defined on the lattice C, in the local basis {\ls) , l^s) , ■ ■ ■ , of V^, can be 
expanded as follows 

d d d 

I*) = E E ■ ■ ■ E ^^1^--^ 1^1) 1^2) • • • Kl) , (64) 

«l=l 42 = 1 ii, = l 

i.e. can be encoded in a tensor T whose d^ elements Ti^i^...i^ G C fully determine the state. 
We refer to the index 1 ^ is ^ d, labeling a local basis for site s, as a physical index. 

The tensor network approach [56] is a powerful way to rewrite the exponentially large 
tensor T in Eq. (64) as a combination of smaller tensors. However, the manipulations 
required to make the products of these tensors easily become too long and cumbersome to 
write explicitly all indices and then it is very convenient to represent them by diagrams. In 
these diagrams each tensor is represented by circles having some outgoing lines and each of 
them represents a tensor's index. For each tensor, we denote its complex conjugate with the 
same circle delimitated by a dashed line instead of a continuos one. A line shared by two 
tensors represents a contraction of a particular index. 

In order to show the power and simplicity of this diagrammatic representation in a 
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T2\2 




(01 



Figure 4. Graphical representation of the identity Tr(/9^^)"'= ~ (Tr/3^° for — 2. 



concrete example, we will graphically prove the identities Tr(p^' 



T2^ 



(Trp 



and 



Tr(pj^2)"° = Trp^° [25], which hold when i? = and the whole chain A = U is in a 
pure state. Let us consider the case of a chain A made by 8 spins and a bipartition such that 
both Ai and A2 contain 4 spins (but the validity of the proof does not rely on these numbers). 
The pure state |\&) is depicted in the top of Fig. 3 where also its Schmidt decomposition is 
shown (the different colors of the circles refer to Ai and A2). The small (purple) circle on 
the line connecting the two big circles stands for the diagonal matrix formed by the Schmidt 



coefficients c^, i.e. |\&) = XlaCalea-*) ® \ea'). We recall that the Schmidt vectors \ea') and 

(2) 

\ea ) provide an orthonormal basis and this ensures that the transformation matrices Ui and 
U2 are isometric [57], a property that graphically can be depicted as in Fig. 3 (top right). 
The density matrix pA = |^)(^| is given by the picture on the bottom (left part) of Fig. 3. 
The definition of given in Eq. (3) can be implemented in the tensor network language by 
interchanging the tensor associated to A2 and its complex conjugate, as again shown in Fig. 3 
(bottom right). 

Now we are ready to show how the above identities for pure states can be easily proved 
by using the graphical representation. We first consider even powers, i.e. Tr(p^^)"'= = 
(Trp^'J,''^)^. In Fig. 4, for sake of graphical simplicity, we restrict to rig = 2 (which trivially 
gives 1 because of the normalization condition TrpA2 = 1) but the key steps hold for any 
even ng. Tr(p^^)^ is obtained by contracting the indices of two copies of as shown in Fig. 
4 (left), and then realizing that the result factorizes into the product of two TrpA2 as in the 
middle of Fig. 4. Using the isometry property in Fig. 3, each Trp^a further reduces to the 
sum of the squared Schmidt coefficients, as it should. Tr(p^^)"° is graphically calculated in 
Fig. 5 for Ho = 3. As above, we first write Tr(pJ)'^ through the proper contractions of three 
copies of p'^ . Then, by rearranging the order of the tensors, this can be put as in the center of 
Fig. 5, which is Trp^^- Again, the last step can be obtained through the isometry property in 
Fig. 3. These two examples make manifest the power of the graphical notation. Finally, it is 
worth stressing, that the same proof could have been done in any other basis. We choose the 



,(2) 
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Figure 5. Graphical representation of the identity Tr(/9^^)"° — Trp^° for no — 3. 



Schmidt one only because in this basis the reduced density matrix is diagonal. 



4.2. The Tree Tensor Network 

Hereafter we focus on the specific tensor network that we have used for the calculation of the 
partially transposed reduced density matrix, i.e. the Tree Tensor Networks (TTN). We only 
briefly recall the basics of the method and for a more detailed explanation, we refer the reader 
to Refs. [55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68]. 

Representing the pure state as in Eq. (64) with a TTN implies that the coefficients 
of the tensor Ti^i^...i^ can be obtained as the result of the contraction of a network made by 
tensors w (built with a smaller number of indices) forming a pattern whose shape reminds a 
tree structure, as graphically shown in Fig. 6 for a lattice of L = 16 sites and tensors w which 
have three indices. The TTN decomposition of Ti^i^...i^ is given by a collection of tensors w, 
hierarchically organized in layers labelled by r, with r = 0, . . . , and (A^ + 1) = log2 L. The 
tensors w have both bond indices and physical indices. All the bond indices are contracted 
while the L physical indices are left non contracted, so that they can be thought as the leaves 
of the tree. In this way, the TTN encodes the d^ complex coefficients Ti-^i^...i^ of the state in 
Eq. (64). 

In the example of Fig. 6, each elementary tensor w has at most one lower index (leg) a 
and two upper indices (legs) /3i and and for this reason is called a binary tree, which is the 
only structure we will consider. At each level r, we allow the indices (3i to run from 1 to x,-. 
Without loss of generality, the tensor w can be chosen to be isometric [57], i.e. = id^. 
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Figure 6. A state encoded by a tree tensor network involving 16 spins. 



with id;^; being the x ^ X identity tensor. Explicitly, this condition reads 



J2 = Saa' , r = 0, . . . , iV - 1 , a, a' = 1, ... , Xr+i • (65) 



/3i,/32=l 

The bottommost tensor of the network (the red one in Fig. 6), i.e. the only tensor at the layer 
A^, is different from the others. Indeed, since it does not possess any lower index, it represents 
a normalized vector, namely 



The properties (65) and (66), which are graphically shown on the right of Fig. 6, guarantee 
that the state |\E') in (64) is normalized to 1, i.e. (^|^) = 1. These considerations can be 
generalized to cases with elementary tensors having more upper and lower legs [58]. The 
natural layered structure of the TTN is emphasized by the arrow on the left of Fig. 6. Each 
layer is an isometric transformation that maps a lattice Cr consisting of sites to a lattice 
Cr+i with Lr+i = Lr/2 sltcs. The physical lattice is Cq, whose length Lq = L must be 
a power of 2. The arrow on the left of Fig. 6 shows the direction along which the coarse- 
graining of the lattice increases. 

In order to describe the ground state of a local Hamiltonian H, the tensors of the 
network should describe the state |\1') that minimizes {^\H\^). This can be achieved using 
the (numerical) algorithm described in Ref. [58] whose computational cost is xtnax^ where 
Xmax = max^lx^} is the largest bond dimension in the network. Recently, the TTN has 
been combined with Monte Carlo sampling obtaining an algorithm whose cost scales with 
Xmax U^]- Each layer of the TTN encoding the ground state of if is a coarse-graining 
transformation that selects those states of Cr which are relevant for the low energy physics 
of H at the scale of Cr+i- When both the Hamiltonian and its ground state are translational 
invariant, we can force the coarse graining transformations to map translational invariant states 
on Cr into translational invariant states on Cr+i by choosing the elementary isometrics of the 
layer r to be all equal. This is represented by the color scheme adopted in Fig. 6, where all 
the isometrics belonging to the same layer have the same color. 




T = N. 



(66) 



/3i,/32 
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Figure 7. Quantities involved in the computation of the pA through a TTN algorithm. Top: A 
is made by two adjacent intervals. Bottom: A is made by two disjoint intervals. 

4.3. Computing the reduced density matrix 

Here we briefly describe how to compute the spectrum of the reduced density matrix pA of a 
subsystem A in a given spin chain [58, 40]. To this aim we should first notice that even when 
the state is translational invariant, by construction the TTN does not enjoy this symmetry. 
As a consequence, there are choices of the subsystem A for which the computation of the 
spectrum of pA is simpler. This happens when A is chosen to be one block (or several blocks 
Ai) of spins of the original lattice which are coarse grained by the TTN to a single spin at 
some level f < N. For example, given the structure in Fig. 6, it is easier to compute the 
spectrum of pA made by two spins when they belong to the same isometry, since they are 
coarse grained to a single spin in one layer. This is due to the fact that, in order to compute 
the spectrum of pA, we do not need to build the full matrix pA, but only a simplified matrix 
Pa which is related to pa by similarity transformation, i.e. 

Pa = Pa ■ (67) 

The matrix VIa corresponds to the contraction of all the isometrics tu involved in this coarse- 
graining process and it transforms each original block Ai into a single spin (see Fig. 7). The 
matrices pa are xT ^ XT^ being Xf the dimensions of the open bonds in pA and m is the 
numbers of disjoint intervals composing A (in this paper we will limit have to m = 1 and 
m = 2). Although pA and pA have different dimensions (with pA smaller than pa), they have 
the same spectrum, because they are related by isometrics. 

Thus, we can focus on the computation of pa- In this process, we can identify three 
groups of tensors: (i) the tensors whose support is fully contained in B; (ii) the tensors whose 
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Figure 8. Quantities involved in the computation of the through a TTN algorithm. Top: 
A is made by two adjacent intervals. Bottom: A is made by two disjoint intervals. 

support is fully contained in A; (iii) the tensors whose support is shared between A and 
its complement. The tensors (i) disappear from the calculation because they are contracted 
through Eq. (65). The tensors (ii) just define ^Ia and therefore they do not occur in the 
computation of pa- Thus, only the tensors of type (iii) are needed for the computation of pa 
which is achieved by contracting these tensors as shown in Fig. 7. Once the simplified tensor 
network encoding pA has been contracted, its diagonalization provides the spectrum of pa- 
The upper bound of the computational cost of the whole algorithm is Xmax- 

Different block configurations (e.g. not optimal choices of Ai and A2 or odd Ai made by 
odd number of spins) can also be obtained from the TTN. Naively, one could think that this 
requires a higher computational cost, but, with some efforts, it is possible to envisage more 
complicated computational strategies to reduce the cost to Xmax^ discussed in Ref. [69] in 
a different context. 

4.4. Computing the partial transpose 

The calculation of the spectrum of the partial transpose of p^ is very similar to the one 
explained above for the spectrum of pa- First we have to choose A as above, namely such 
that each Ai becomes a single spin at level f. Then, we can identify three groups of tensors : 
(I) the tensors fully contained into B; (II) the tensors fully contained into either Ai or A2 (not 
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both); (III) the tensors connecting Ai and A2 or Ai and B. Notice that, while (I) is the same 
as (i) above, the group (II) is different from (ii) and, consequently, also (III) is not (iii). Let us 
focus on computing 

pl- = n'AP^l'^n'\, (68) 

where Q'a is an isometric tensor (i.e. Q'a^'a — id^^t")- For the same reasons discussed in the 

~(2) 

subsection 4.3, only the tensor of type (III) enter in the definition of . 

~(2) 

When A is made by two disjoint blocks (see Fig. 8, bottom panel), the matrix p\' is 
given by = p J . For two adjacent blocks (see Fig. 8, top panel) this relation does not hold 
because pa has only one incoming and one outgoing bond and therefore we cannot distinguish 
between Ai and A2, which is necessary in order to perform the partial transposition. In this 
case, indeed we have ^I'a 7^ and the former reads 

n'A = nA,0nA,. (69) 

Comparing the top panel of Fig. 8 (center) to the top panel of Fig. 7 (center), we realize that 
the yellow isometry connecting Ai to A2 does not occur in Q'a, in contrast with f2^. Indeed, 

~(2) 

the yellow isometry appears in p^ . We finally mention that the computational cost of all the 
algorithm to determine the spectrum of pj is upper-bounded by Xmax- 



5. Numerical results for the critical Ising chain in a transverse magnetic field 

In this section we report numerical results for the transverse field Ising chain obtained by 
means of the tree tensor network as explained in the previous section. We consider chains of 
finite length L of the form L = 2^^ and the maximum size studied is L = 512. The Xmax 
for this simulation has been fixed to 128, guaranteeing a relative precision on the ground state 
energy of about 10^*^. We recall that with the TTN method, using a binary tree, we can quickly 
access only subsystems of size £ = 2™ with m < M integer, as it should be clear from the 
previous section. 

As a first calculation we considered a bipartite chain (i.e. 5—7-0) and we checked 
that the entanglement negativity reproduces the Renyi entropy S^^^'^\ as it should. We do 
not report explicit plots of these results, but they have been important numerical checks for 
the numerics. The results for tripartite systems are reported and discussed in the following 
subsections. 



5.1. Two adjacent intervals. 

We first consider the case of two adjacent intervals both of equal length ^ in a periodic chain 
of total length L so that all the results depend on the single parameter z = i / L E [0, 1/2]. In 
terms of z and fixing c = 1/2, the CFT predictions in Eq. (29) can be written as 

(L/tt sin(7r2))-("^/2-2/"^)/6(L/7r sin(27r2))-("^/2+i/ne)/i2 ^ 
Tr(pJ)" = dnx{ (70) 
((L/7r)3 sin^nz) sin(27rz))-("°-i/"°)/24, 
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Figure 9. Entanglement for two adjacent intervals of equal length £ < L/2 in a periodic 
chain of length L. The quantity r„(z) in Eq. (72) as function of z = £/L compared with the 
parameter free CFT prediction for n = 3 (left) and n ~ A (right). 



while for the logarithmic negativity we have from Eq. (30) 

L 



^ = - hi 



7r 



+ cnst. (71) 



Following Ref. [25], we can construct quantities in which the dependence on the non- 
universal parameters dn and also the universal dependence on L cancel. To this aim, it is 
enough to divide Tr(p^^)" by the value it assumes at a given fixed ^, e.g. ^ = L/A, i.e. by 
considering the quantities 

r^{z) = In ^'^^f 2" ' (^2) 
whose parameter free CFT predictions for n even and odd are 

'"''='"''^-"'> -U^^^) 

'^no = 777 f — - '^o) ln(2 sin^(7rz) sin(27r2)). (73) 

The numerical results for these quantities are shown in Fig. 9 for n = 3 and n = A. The 
agreement between the numerical data and the CFT predictions is perfect for all considered 
values of L, showing that finite size corrections are very small for these quantities. Notice that 
for 2; = 1/2 we have a bipartite system (i.e. i? -> 0) and the equations in (73) obviously do 
not work since the data are described by Eqs. (26) and (27), reflecting the fact that the limit 
z — i- 1/2 is approached in a non-uniform way (i.e. the limits z ^ 1/2 and N ^ 00 Ao not 
commute as obvious). 

For the logarithmic negativity, we can analogously define the subtracted quantity 

e{z) = S{i, L) - £{L/A, L) = \ ln[tan(7rz)] , (74) 

8 

and again the r.h.s. is a parameter free CFT prediction. In Fig. 10, this prediction is compared 
with the numerical data and the agreement is extremely good except at 2 = 1/2 where the 
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Figure 10. Entanglement negativity for two adjacent intervals of equal length £ < L/2 in a 
periodic chain of length L: Subtracted negativity e{z) in Eq. (74) compared with the parameter 
free CFT prediction. 



numerical data are described by the bipartite formula (28). Notice that finite size scaling 
corrections are even smaller than those for the quantities r„(2;) in Fig. 9. 

5.2. Two disjoint intervals 

In this section we study the most interesting and difficult situation of two disjoint intervals for 
which an accurate numerical study of the negativity has been already performed by means of 
density matrix renormalization group in Ref. [15], but before the systematic CFT derivation 
in Refs. [24, 25]. Here we first consider the traces Tr(p^^)" and Tr(pA)". Indeed, although 
standard Renyi entropies given by Trp^ have been already studied in Refs. [30, 38, 39], they 
show large corrections to the scaling which is worth recalhng before embarking in the analysis 
ofTr(pJ)". 

In order to determine numerically the function we calculate for several finite 

chains the quantity 

F':\x) = (1 - , (75) 

which in the scaling limit should converge to the CFT prediction Tnix). In a finite system 
of length L, the four-point ratio x must be rewritten by replacing all distances by the 
corresponding chordal lengths. For two intervals of the same length i at distance r this reads 

v. (76) 



^sin(7r(£ + r)/L) 

The numerical data for the function F^^{x) are reported in Fig. 11 for n = 3,4,5 
as function of x for various values of i (i.e. different values of L according to Eq. 
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Figure 11. Numerical results for Fl^^{x) as function of x for different values of ^ for n = 3 
(top left), n — 4 (top right), and n = 5 (bottom). The data are extrapolated to £ — > oo by 
means of Eq. (78). The extrapolated data (topmost set of data) are in excellent agreement with 
the CFT prediction (continuous line). 



(76)). It is evident that strong scaling corrections affect the data, as known from previous 
analyses [38, 39]. It has been argued on the basis of the general CFT arguments [72], and 
shown explicitly in few examples [71, 73, 74] both analytically and numerically, that the 
entanglement entropies (hence also the function F^^(x)), display 'unusual' corrections to the 
scaling which, at the leading order, can be effectively taken into account by the scaling ansatz 

F'fix) = J^nix) + r'/^fn{x) + .... (77) 

However, the corrections to the scaling in Fig. 11 cannot be captured by this simple 
ansatz because subleading corrections to the scaling become more and more important with 
increasing the value of n. Indeed, corrections of the form for any integer m are know 

to be present [73, 39, 40]. Thus the most general finite-^ ansatz is 

J^"W = ^„W + 4^ + ^|M + ^.... (78) 

The effect of the subleading corrections is enhanced by the fact the amplitude functions fn^ (x) 
have different signs determining a non-monotonic behavior in £ (in particular, we have that 
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Figure 12. Numerical results for G^* {x) as function of y for different values of £ for n — 3 
(top left), 71 = 4 (top right), and n = 5 (bottom). The data are extrapolated to £ — > oo by 
means of Eq. (80). The extrapolated data (topmost set of data) are in excellent agreement with 
the CFT prediction (continuous line). 



fn'^ is always negative, while fn^ and fn^ are always positive, as discussed in Ref. [40]). In 
order to have an accurate extrapolation to £ — > oo, for any n we consider all the corrections 
above up to order and we get the extrapolations reported in Fig. 11. The error bars 

are estimated by studying the stability of the extrapolation with respect to the number of sizes 
i included in the fit. The overall agreement of the extrapolated points with the CFT prediction 
is excellent for all values of x and for the three considered values of n, reproducing the results 
inRefs. [38,39]. 

After having summarized the corrections to the scaling for the entanglement entropies 
we can turn to the integer powers of the partial transpose in which we are interested here. In 
analogy with Eq. (75) we can define the lattice ratio 

G':\y) = P^^il - yr-'/-^/'' , (79) 

that in the limit £ — > oo is expected to converge to the CFT scaling function Qniv) given by 
Eq. (53). In the case at hand, the numerical value of y is given by the same expression in Eq. 
(76) for X. The numerical data for Gl^^{y) are reported in Fig. 12 forn = 3, 4, 5 as function 
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Figure 13. Numerical results for i?Jf'(a;) as function of y for different values of £ for n — 3 
(top left), n — 4 (top right), and n = 5 (bottom). The data are extrapolated to £ — > oo by 
means of Eq. (82). The extrapolated data (topmost set of data) are in excellent agreement with 
the CFT prediction (continuous line). 



of y for various values of i. Also in this case, large scaling corrections are present. These are 
expected to be of the same form as for i^^'^*(a;), i.e. described by the ansatz 

[y) = Gniy) + + + .... (80) 

We repeat exactly the same analysis as for F^'^*(a;) to extrapolate the data to £ — )■ cxd. As above, 
we find that giP is always negative, while Qn^ and gn^ are positive. Using this observation, we 
extrapolate to £ — )■ oo, obtaining the results (with error bars) reported in Fig. 12. These points 
agree very well with the CFT prediction in Eq. (53) for the three values of n considered. 

It is evident from the comparison of Figs. 1 1 and 1 2 that the functions Gjf* (y) and F^^* (y) 
are very close to each other and the logarithmic negativity is instead only determined by the 
small differences between the two. Thus, as already discussed in Sec. 3, a practical way to 
have more accurate tests of the CFT predictions which are sensitive to the small differences 
between the two universal functions Gniv) and -F„(y) is to consider the ratio Rniy) between 



Entanglement negativity in the critical Ising chain 



26 



0.03 
0.025 

0.02 
0.015 

0.01 
0.005 






1 = 


8 




1 = 


16 


+ 


1 = 


32 




1 = 


64 





1 = 


128 







0.1 



0.2 0.3 

Y 



0.4 0.5 



Figure 14. Logarithmic negativity for two intervals of equal length t at distance r as function 
of the four point ratio y. 



the two and the analogous lattice quantity 



Rn\y) 



(81) 



which in the limit £ — )■ oo converges to the CFT prediction in Eq. (54). The numerical data 
for Rn^{y) are reported in Fig. 13 for n = 3, 4, 5 as function of y for different values of i. 
Once again, large scaling corrections are present and there are no accidental cancellations in 
the ratio, so that they are again expected to be of the same form as for F^^{x), i.e. described 
by the ansatz 



rn\y) ^ r)t'{y) 



(2), 



+ 



^(3) 
' n 



(82) 



£l/n £2/n ' £3/n 

We repeat again the same analysis as for F^^{x) to extrapolate the data to £ — )• oo and the 
results (with error bars) are reported in Fig. 13. Unlike fn^ (s) 's and gn^ (?/) 's, in this case the 

(i) / \ (?) 

signs of r„ (y)'s are not defined (indeed rn s can be written as complicated combinations of 
fn^'s and Qn^'s). For this reason, the error bars in Fig. 13 are larger than the ones in Fig. 1 1 
and in Fig. 12. It is evident that the extrapolated points in Fig. 13 agree very well with the 
CFT prediction for the three considered values of n. It is very remarkable that the numerical 
calculations are accurate enough to detect the small differences of these ratios from 1 (at least 
for n = 3 and n = 4, while for n = 5 the estimated error is too large to distinguish the 
extrapolation from one). 

Finally we turn to the study of the logarithmic negativity £. The numerical data as a 
function of y are reported in Fig. 14 for several values of £. In the figure all data collapse on 
a single curve, with some tiny corrections to the scaling for the smaller values of i, which 
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however are much smaller than the scaling corrections for the ratios Rniy), as evident from 
a comparison between Fig. 14 and Fig. 13. This is similar to what has been observed for 
the free boson in Ref. [25] and recalls the similar effect for the entanglement entropies both 
for one and two intervals where the unusual corrections to the scaling are suppressed for the 
von Neumann entropy for periodic boundary conditions [71, 73]. It would be then interesting 
to study the negativity for Ising chains with open boundary conditions to check whether the 
unusual corrections are enhanced (analogously to Friedel oscillations [75, 76]). On the more 
theoretical side, we have not been able to find the analytic continuation of Rn^ (y) to ne — )• 1 
which would allow a strict check of the numerical data in Fig. 14. However, the data for 
y <^ 1 are very close to zero in agreement with the general result [25] that the negativity 
should fall off faster than any power with the separation of the intervals. 

6. Conclusions 

In this manuscript we presented a systematic study of the entanglement negativity and of 
the traces of integer powers of the partial transpose of the reduced density matrix in the 
critical Ising chain. In order to have an accurate numerical evaluation of these quantities we 
adapted the tree tensor network technique to the calculation of the eigenvalues of the partially 
transposed reduced density matrix. For two adjacent intervals we found perfect agreement 
with the model independent CFT predictions in Refs. [24, 25]. For two disjoint intervals 
we first derived a CFT prediction for Tr(p^^^^^^)" which later has been compared with the 
numerical data taking into account the finite size corrections induced by the finite length of 
the blocks. Unfortunately, in order to determine the logarithmic negativity, it remains a hard 
open problem to find the analytic continuation to Ue I of these traces. This reflects the 
similar problem found for the standard entanglement entropies, where explicit formulas for 
Trp^juA2 ^^^^ been obtained for some models [29, 30], but they have not been analytically 
continued to n — )■ 1, preventing us to write down a close formula for the von Neumann 
entropy. 
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